**
This R markdown document (.Rmd) contains analysis, visualisation and statistical testing relating to research paper “Reticular adhesions mediate cell-matrix attachment during mitosis” by Lock et al (NCB, 2018).
This code uses the R markdown framework and enables production of statistical and graphical outputs equivalent to those presented in the paper. Some purely aesthetic differences exist to final paper figures given subsequent editing (adobe illustrator) for figure embedding. All R packages necessary to run this code will be installed automatically upon running this code.
This code calls the source data file named “Supplementary Table 1 Statistics Source Data.xlsx” and selects individual data sheets pertinent to each analysis. These sheets are named according the Figure Panel that they underpin. The user must place this source data file in a location of their choice and modify the “data_input_path” variable to reflect this location. Graphical outputs are saved according to the relevant Figure Panel names (with some additional descriptive information) to a location defined by the variable “panel_output_path”. By default this is in a folder named “Graphical_Outputs” within the defined “data_input_path”. Note that some figure panels were generated directly within Excel and these can be found within the source data file.
In addition to specific .svg graphical outputs, this R markdown code will produce an HTML notebook that is useful to collectively preview and efficiently share outputs and their underlying code. This output HTML notebook file matches the name of this source .Rmd file, i.e. “Statistical_Analysis_Code_Lock_et_al_NCB_2018.nb.html”.
**
rr # Set your path to input data data_input_path <- /Users/johnlock/Dropbox/A-CMAC Manchester collaboration/NCB submission documents/Resubmission 2/Source_Data
# name of the source data file which you have placed in the above location data_file <- Table 1 Statistics Source Data.xlsx
rr # Define panel output path panel_output_path <- paste(data_input_path, _Outputs/, sep = /)
rr # Import individual data files Mass.Spec.Data <- read_xlsx(paste(data_input_path, data_file, sep = /), sheet = 1A, col_names = TRUE) b5.int.vs.VN.conc <- read_xlsx(paste(data_input_path, data_file, sep = /), sheet = 1G, col_names = TRUE)
rr #define panel_file_name Integrin_Spectral_Counts_File_Name <- 1A Integrin Spectral Counts Boxplot
boxplot(Mean.Spectral.Counts ~ Integrin.Subunit, data = Mass.Spec.Data, notch = FALSE, varwidth = FALSE, las=1, ylab = Spectral Counts, col = , ylim = c(0, 50), outline = FALSE, boxlty = 0, whisklty = 0, staplelty = 0, medlwd = 0.5, boxfill = NA) beeswarm(Mean.Spectral.Counts ~ Integrin.Subunit, data = Mass.Spec.Data, method = ‘center’, add = T, col = , pch = 19, cex = 0.7) r svglite(paste(panel_output_path, Integrin_Spectral_Counts_File_Name, .svg, sep = \), width=8, height=6) boxplot(Mean.Spectral.Counts ~ Integrin.Subunit, data = Mass.Spec.Data, notch = FALSE, varwidth = FALSE, las=1, ylab = Spectral Counts, col = , ylim = c(0, 50), outline = FALSE, boxlty = 0, whisklty = 0, staplelty = 0, medlwd = 0.5, boxfill = NA) r beeswarm(Mean.Spectral.Counts ~ Integrin.Subunit, data = Mass.Spec.Data, method = ‘center’, add = T, col = , pch = 19, cex = 0.7) dev.off()
quartz_off_screen
2
rr AVvsB1 <- subset(Mass.Spec.Data, Integrin.Subunit == | Integrin.Subunit == 1) T_test_AVvsB1 <- t.test(Mean.Spectral.Counts ~ Integrin.Subunit, data = AVvsB1, alternative = .sided, paired = FALSE)$p.value T_test_AVvsB1
[1] 0.005444416
rr B5vsB1 <- subset(Mass.Spec.Data, Integrin.Subunit == 5 | Integrin.Subunit == 1) T_test_B5vsB1 <- t.test(Mean.Spectral.Counts ~ Integrin.Subunit, data = B5vsB1, alternative = .sided, paired = FALSE)$p.value T_test_B5vsB1
[1] 0.0001080229
rr #define panel_file_name b5_Intensity_File_Name <- 1G Integrin b5 Intensity
names(b5.int.vs.VN.conc) <- c(_Pos_vs_Neg, , .integrin.b5.Intensity) b5.int.vs.VN.conc\(Talin_Pos_vs_Neg <- recode(b5.int.vs.VN.conc\)Talin_Pos_vs_Neg, Atypical = _Neg) b5.int.vs.VN.conc\(Talin_Pos_vs_Neg <- recode(b5.int.vs.VN.conc\)Talin_Pos_vs_Neg, typical = _Pos) b5.int.vs.VN.conc\(Group <- paste(b5.int.vs.VN.conc\)Talin_Pos_vs_Neg, b5.int.vs.VN.conc\(Condition, sep = \_\) b5.int.vs.VN.conc\)Group <- factor(b5.int.vs.VN.conc$Group, levels = c(_Pos_1VN, _Neg_1VN, _Pos_3VN, _Neg_3VN, _Pos_10VN, _Neg_10VN)) ggplot(data = b5.int.vs.VN.conc, aes(x = Group, y = Mean.integrin.b5.Intensity)) + geom_boxplot(aes(fill = Talin_Pos_vs_Neg, ymin=..lower.., ymax=..upper..), outlier.shape=NA, notch = TRUE) + scale_fill_manual(values = c(_Pos = #4F80BC, _Neg = #C0504F)) + theme_classic() + coord_cartesian(ylim = c(500, 1300)) + theme(axis.text.x = element_text(angle = 45, hjust = 1)) svglite(paste(panel_output_path, b5_Intensity_File_Name, .svg, sep = \), width=5, height=8) r ggplot(data = b5.int.vs.VN.conc, aes(x = Group, y = Mean.integrin.b5.Intensity)) + geom_boxplot(aes(fill = Talin_Pos_vs_Neg, ymin=..lower.., ymax=..upper..), outlier.shape=NA, notch = TRUE) + scale_fill_manual(values = c(_Pos = #4F80BC, _Neg = #C0504F)) + theme_classic() + coord_cartesian(ylim = c(500, 1300)) + theme(axis.text.x = element_text(angle = 45, hjust = 1)) dev.off()
quartz_off_screen
2
rr TalinPos = subset(b5.int.vs.VN.conc, b5.int.vs.VN.conc\(Talin_Pos_vs_Neg == \Talin_Pos\) TalinPos_VN1_vs_VN3 <- subset(b5.int.vs.VN.conc, b5.int.vs.VN.conc\)Talin_Pos_vs_Neg == _Pos & b5.int.vs.VN.conc\(Condition == \1VN\ | b5.int.vs.VN.conc\)Talin_Pos_vs_Neg == _Pos & b5.int.vs.VN.conc\(Condition == \3VN\ ) TalinPos_VN3_vs_VN10 <- subset(b5.int.vs.VN.conc, b5.int.vs.VN.conc\)Talin_Pos_vs_Neg == _Pos & b5.int.vs.VN.conc\(Condition == \3VN\ | b5.int.vs.VN.conc\)Talin_Pos_vs_Neg == _Pos & b5.int.vs.VN.conc\(Condition == \10VN\ ) Mu_Test_TalinPos_VN1_vs_VN3 = wilcox.test(Mean.integrin.b5.Intensity ~ Condition, data = TalinPos_VN1_vs_VN3)\)p.value Mu_Test_TalinPos_VN3_vs_VN10 = wilcox.test(Mean.integrin.b5.Intensity ~ Condition, data = TalinPos_VN3_vs_VN10)\(p.value TalinNeg_VN1_vs_VN3 <- subset(b5.int.vs.VN.conc, b5.int.vs.VN.conc\)Talin_Pos_vs_Neg == _Neg & b5.int.vs.VN.conc\(Condition == \1VN\ | b5.int.vs.VN.conc\)Talin_Pos_vs_Neg == _Neg & b5.int.vs.VN.conc\(Condition == \3VN\ ) TalinNeg_VN3_vs_VN10 <- subset(b5.int.vs.VN.conc, b5.int.vs.VN.conc\)Talin_Pos_vs_Neg == _Neg & b5.int.vs.VN.conc\(Condition == \3VN\ | b5.int.vs.VN.conc\)Talin_Pos_vs_Neg == _Neg & b5.int.vs.VN.conc\(Condition == \10VN\ ) Mu_Test_TalinNeg_VN1_vs_VN3 = wilcox.test(Mean.integrin.b5.Intensity ~ Condition, data = TalinNeg_VN1_vs_VN3)\)p.value Mu_Test_TalinNeg_VN3_vs_VN10 = wilcox.test(Mean.integrin.b5.Intensity ~ Condition, data = TalinNeg_VN3_vs_VN10)$p.value
rr # Set your path to input data data_input_path <- /Users/johnlock/Dropbox/A-CMAC Manchester collaboration/NCB submission documents/Resubmission 2/Source_Data
data_file <- Table 1 Statistics Source Data.xlsx
# Import individual data files FRAP.data <- read_xlsx(paste(data_input_path, data_file, sep = /), sheet = 2S_T, col_names = TRUE) STORM.Data <- read_xlsx(paste(data_input_path, data_file, sep = /), sheet = 2V_W, col_names = TRUE)
rr #define panel_file_name FRAP_Curve_File_Name <- 2S b5 FRAP recovery
names(FRAP.data) <- c(_b5_FRAP_recovery_relative_to_prebleach, , _ROI, _vs_Focal) FRAP.data\(Reticular_vs_Focal <- recode(FRAP.data\)Reticular_vs_Focal, atypical = ) FRAP.data\(Reticular_vs_Focal <- recode(FRAP.data\)Reticular_vs_Focal, typical = ) ggplot(data=FRAP.data, aes(x=Time, y=Standardised_b5_FRAP_recovery_relative_to_prebleach, colour = Reticular_vs_Focal), alpha = 0.1) + stat_summary(fun.data =_cl_normal, alpha = 0.5) + theme_classic() + scale_colour_manual(values = c( = #BF4C49, = #4B7FBB)) + geom_smooth(method = , span = 0.2, size = 1, fill = #3A3A3A) svglite(paste(panel_output_path, FRAP_Curve_File_Name, .svg, sep = \), width=6, height=4) r ggplot(data=FRAP.data, aes(x=Time, y=Standardised_b5_FRAP_recovery_relative_to_prebleach, colour = Reticular_vs_Focal), alpha = 0.1) + stat_summary(fun.data =_cl_normal, alpha = 0.5) + theme_classic() + scale_colour_manual(values = c( = #BF4C49, = #4B7FBB)) + geom_smooth(method = , span = 0.2, size = 1, fill = #3A3A3A) dev.off()
quartz_off_screen
2
rr Reticular.FRAP.data <- subset(FRAP.data, Reticular_vs_Focal == ) Focal.FRAP.data <- subset(FRAP.data, Reticular_vs_Focal == ) Reticular.FRAP.data.LOESS <- loess(Standardised_b5_FRAP_recovery_relative_to_prebleach ~ Time, data = Reticular.FRAP.data, span = 0.2) Focal.FRAP.data.LOESS <- loess(Standardised_b5_FRAP_recovery_relative_to_prebleach ~ Time, data = Focal.FRAP.data, span = 0.2) Reticular.FRAP.data.LOESS.singular <- unique(Reticular.FRAP.data.LOESS\(fitted) Focal.FRAP.data.LOESS.singular <- unique(Focal.FRAP.data.LOESS\)fitted) KS.Reticular.vs.Focal.FRAP <- ks.test(Reticular.FRAP.data.LOESS.singular, Focal.FRAP.data.LOESS.singular, alternative = .sided) print(KS.Reticular.vs.Focal.FRAP$p.value)
[1] 0.002107921
rr ## Determine CMAC.type (Reticular vs Focal depending on presence of A in \(File (column)) CMAC.type = c() for (cluster in 1:length(STORM.Data\)File) ) { CMAC.type[cluster] <- if(grepl(, STORM.Data\(File[cluster]) == TRUE) { \Reticular\
} else if(grepl(\T\, STORM.Data\)File[cluster]) == TRUE) {
} else if(grepl(, STORM.Data\(File[cluster]) == TRUE) { \Non-Retraction\
} else { \Retraction\
} } STORM.Data\)CMAC.type <- cbind(as.character(CMAC.type)) ## Unique nanocluster ID generation = merge of Folder, File and Cluster # STORM.Data\(Unique.ncID <- as.factor(paste(STORM.Data\)Folder, STORM.Data\(File, STORM.Data\)Cluster)) STORM.Data\(Unique.CMACID <- as.factor(paste(STORM.Data\)Folder, STORM.Data$File)) STORM.Data <- as.data.frame(STORM.Data)
rr ## NND analysis of Incite3 data NNND <- c() for (i in 1:length(unique(STORM.Data\(Unique.CMACID)) ) { temp = STORM.Data[STORM.Data\)Unique.CMACID == unique(STORM.Data$Unique.CMACID)[i], ] NNND = append(NNND, nndist(temp[,7], temp[,8], k=1)) }
NAs introduced by coercionNAs introduced by coercion
rr STORM.Data$NNND <- cbind(NNND)
rr NNND.boxplot.filename <- 2V nearest neighbour boxplot
NNND.A = subset(STORM.Data\(NNND, CMAC.type == \Reticular\) NNND.T = subset(STORM.Data\)NNND, CMAC.type == ) svglite(paste(panel_output_path, NNND.boxplot.filename, .svg, sep = \), width=6, height=4) boxplot(NNND ~ CMAC.type, subset(STORM.Data, CMAC.type == | CMAC.type == ), notch = TRUE, outline = FALSE, main = Neighbour Distance of Nanocluster by CMAC Type, xlab = type, ylab = Neighbour Distance (nm), plot =
) dev.off()
null device
1
rr boxplot(NNND ~ CMAC.type, subset(STORM.Data, CMAC.type == | CMAC.type == ), notch = TRUE, outline = FALSE, main = Neighbour Distance of Nanocluster by CMAC Type, xlab = type, ylab = Neighbour Distance (nm), plot =
)
rr Molecule.number.boxplot.filename <- 2W molecular number boxplot
Molecules.A = as.numeric(subset(STORM.Data\(Molecules, CMAC.type == \Reticular\)) Molecules.T = as.numeric(subset(STORM.Data\)Molecules, CMAC.type == )) svglite(paste(panel_output_path, Molecule.number.boxplot.filename, .svg, sep = \), width=6, height=4) boxplot(as.numeric(Molecules) ~ CMAC.type, subset(STORM.Data, CMAC.type == | CMAC.type == ), notch = TRUE, outline = FALSE, main = of Integrin avb5 Molecules by CMAC Type, xlab = type, ylab = of Integrin avb5 Molecules, plot =
) dev.off()
null device
1
rr boxplot(as.numeric(Molecules) ~ CMAC.type, subset(STORM.Data, CMAC.type == | CMAC.type == ), notch = TRUE, outline = FALSE, main = of Integrin avb5 Molecules by CMAC Type, xlab = type, ylab = of Integrin avb5 Molecules, plot =
)
rr # Set your path to input data data_input_path <- /Users/johnlock/Dropbox/A-CMAC Manchester collaboration/NCB submission documents/Resubmission 2/Source_Data
data_file <- Table 1 Statistics Source Data.xlsx
# Import individual data files CytoD2h_vs_control <- read_xlsx(paste(data_input_path, data_file, sep = /), sheet = 3C, col_names = TRUE) CytoD2h_vs_control <- CytoD2h_vs_control[,1:3] Adhesion.Assay.Data <- read_xlsx(paste(data_input_path, data_file, sep = /), sheet = 3E, col_names = TRUE) Adhesion.Assay.Data <- Adhesion.Assay.Data[,1:4] Talin.KD.Data <- read_xlsx(paste(data_input_path, data_file, sep = /), sheet = 3H_I, col_names = TRUE)
rr oldpar <- par() par(mar = c(10, 2.5, 1, 1)) #define panel_file_name Adhesion_Assay_panel_file_name <- 3E Adhesion Assay Boxplots
Adhesion.Assay.Data\(Condition <- factor(Adhesion.Assay.Data\)Condition, unique(Adhesion.Assay.Data\(Condition)) max_adhesion <- max(Adhesion.Assay.Data\)NumberObjects) Adhesion.Assay.Data\(Normalised_Adhesion <- Adhesion.Assay.Data\)NumberObjects/max_adhesion svglite(paste(panel_output_path, Adhesion_Assay_panel_file_name, .svg, sep = \), width=7.7, height=6.5) boxplot(Normalised_Adhesion ~ Condition, data = Adhesion.Assay.Data, notch = FALSE, varwidth = FALSE, outline = F, las=2, ylab = Cell Number, col = c(1wt -Cyto = #CB2327, 1wt +Cyto +RGD = #CB2327, 1b5 -Cyto = #4775A3, 1b5 -Cyto +RAD = #4775A3, 1b5 -Cyto +RGD = #9EC0E5, 1b5 +Cyto = #4775A3, 1b5 +Cyto +RAD = #4775A3, 1b5 +Cyto +RGD = #9EC0E5), ylim = c(0,1.2)) beeswarm(Normalised_Adhesion ~ Condition, data = Adhesion.Assay.Data, method = ‘center’, add = T, col = , pch = 19, cex = 0.7) dev.off()
quartz_off_screen
2
rr boxplot(Normalised_Adhesion ~ Condition, data = Adhesion.Assay.Data, notch = FALSE, varwidth = FALSE, outline = F, las=2, ylab = Cell Number, col = c(1wt -Cyto = #CB2327, 1wt +Cyto +RGD = #CB2327, 1b5 -Cyto = #4775A3, 1b5 -Cyto +RAD = #4775A3, 1b5 -Cyto +RGD = #9EC0E5, 1b5 +Cyto = #4775A3, 1b5 +Cyto +RAD = #4775A3, 1b5 +Cyto +RGD = #9EC0E5), ylim = c(0,1.2)) beeswarm(Normalised_Adhesion ~ Condition, data = Adhesion.Assay.Data, method = ‘center’, add = T, col = , pch = 19, cex = 0.7)
rr p_values <- c() for (i in 1:length(unique(Adhesion.Assay.Data\(Condition))){ for (j in 1:length(unique(Adhesion.Assay.Data\)Condition))){ Condition_1 = subset(Adhesion.Assay.Data, Adhesion.Assay.Data\(Condition == unique(Adhesion.Assay.Data\)Condition)[i]) Condition_1\(Condition_Number = as.factor(i) Condition_2 = subset(Adhesion.Assay.Data, Adhesion.Assay.Data\)Condition == unique(Adhesion.Assay.Data\(Condition)[j]) Condition_2\)Condition_Number = as.factor(j+100) Condition_Pair = rbind(Condition_1, Condition_2) Condition_Pair_ttest = t.test(Normalised_Adhesion ~ Condition_Number, data = Condition_Pair, paired = FALSE)\(p.value p_values = rbind(p_values, Condition_Pair_ttest) } } p_values <- as.vector(p_values) p_value_matrix <- matrix(p_values,nrow = length(unique(Adhesion.Assay.Data\)Condition)),ncol = length(unique(Adhesion.Assay.Data\(Condition))) columns <- as.vector(unique(Adhesion.Assay.Data\)Condition)) rows <- as.vector(unique(Adhesion.Assay.Data\(Condition)) p_value_dataframe <- as.data.frame(p_value_matrix, row.names = rows) names(p_value_dataframe) <- columns # print(p_value_dataframe) # Introduce Holm (also called Holm-Bonferroni) p-value correction (used because it gives same stringency against false positives (type 1 errors) with lower probability of false negatives (type 2 errors))) corrected_p_values <- p.adjust(p_values, method = \holm\, n = length(p_values)) corrected_p_value_matrix1 <- matrix(corrected_p_values,nrow = length(unique(Adhesion.Assay.Data\)Condition)),ncol = length(unique(Adhesion.Assay.Data\(Condition))) columns <- as.vector(unique(Adhesion.Assay.Data\)Condition)) rows <- as.vector(unique(Adhesion.Assay.Data$Condition)) corrected_p_value_dataframe <- as.data.frame(corrected_p_value_matrix1, row.names = rows) names(corrected_p_value_dataframe) <- columns print(corrected_p_value_dataframe)
rr heatmap.2(log(corrected_p_value_matrix1), dendrogram = ‘none’, Rowv = FALSE, Colv = FALSE, labRow = rows, labCol = columns, srtCol = 45, cexRow = 1, cexCol = 1, margins = c(7, 9), symm = TRUE, revC = FALSE, breaks = c(log(1e-300), log(1e-10), log(1e-6), log(0.001), log(0.01), log(0.05), log(1)), rowsep = 1:nrow(corrected_p_value_matrix1), colsep = 1:ncol(corrected_p_value_matrix1), sepcolor = ‘darkgrey’, sepwidth = c(0.02,0.02), trace = ‘none’, col = c(, , , , , ), denscol = NULL, keysize = 1.5, key.title = NA)
rr U_values <- c() for (i in 1:length(unique(Adhesion.Assay.Data\(Condition))){ for (j in 1:length(unique(Adhesion.Assay.Data\)Condition))){ Condition_1 = subset(Adhesion.Assay.Data, Adhesion.Assay.Data\(Condition == unique(Adhesion.Assay.Data\)Condition)[i]) Condition_1\(Condition_Number = as.factor(i) Condition_2 = subset(Adhesion.Assay.Data, Adhesion.Assay.Data\)Condition == unique(Adhesion.Assay.Data\(Condition)[j]) Condition_2\)Condition_Number = as.factor(j+100) Condition_Pair = rbind(Condition_1, Condition_2) Condition_Pair_Utest = wilcox.test(Normalised_Adhesion ~ Condition_Number, data = Condition_Pair)$p.value U_values = rbind(p_values, Condition_Pair_Utest) } }
cannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with tiescannot compute exact p-value with ties
rr U_values <- as.vector(U_values) U_value_matrix <- matrix(U_values, nrow = length(unique(Adhesion.Assay.Data\(Condition)),ncol = length(unique(Adhesion.Assay.Data\)Condition))) columns <- as.vector(unique(Adhesion.Assay.Data\(Condition)) rows <- as.vector(unique(Adhesion.Assay.Data\)Condition)) U_value_dataframe <- as.data.frame(U_value_matrix, row.names = rows) names(U_value_dataframe) <- columns print(U_value_dataframe)
rr # beeswarm(Total.Talin.per.Cell ~ siRNA..Oligo, data = Talin.KD.Data, method = ‘center’, col = , pch = 19, cex = 1, ylim = c(0,1)) # # beeswarm(Cell.Area ~ siRNA..Oligo, data = Talin.KD.Data, method = ‘center’, col = , pch = 19, cex = 1, ylim = c(0,1)) # # beeswarm(Mean.Clustered.b5.per.Cell ~ siRNA..Oligo, data = Talin.KD.Data, method = ‘center’, col = , pch = 19, cex = 1)#, ylim = c(0,1)) plot(log2(Talin.KD.Data\(Total.Talin.per.Cell), Talin.KD.Data\)Cell.Area, ylim = c(0,1), pch = c(0, 0, 0, 0, 2, 21), col = c(, , , , , )) abline(lm(Talin.KD.Data\(Cell.Area ~ log2(Talin.KD.Data\)Total.Talin.per.Cell)), col = ) r lm.cell.area_vs_talin.per.cell <- lm(Talin.KD.Data\(Cell.Area ~ log2(Talin.KD.Data\)Total.Talin.per.Cell)) TalinKD_Cell_Area_panel_file_name <- 3H TalinKD_vs_Cell_Area
svglite(paste(panel_output_path, TalinKD_Cell_Area_panel_file_name, .svg, sep = \), width=3, height=3) plot(log2(Talin.KD.Data\(Total.Talin.per.Cell), Talin.KD.Data\)Cell.Area, ylim = c(0,1), pch = c(0, 0, 0, 0, 2, 21), col = c(, , , , , ), cex.axis = 0.5, cex.label = 0.5) r abline(lm(Talin.KD.Data\(Cell.Area ~ log2(Talin.KD.Data\)Total.Talin.per.Cell)), col = ) dev.off()
quartz_off_screen
2
rr plot(Talin.KD.Data\(Total.Talin.per.Cell, log2(Talin.KD.Data\)Mean.Clustered.b5.per.Cell), xlim = c(0,1), ylim = c(0,4), pch = c(0, 0, 0, 0, 2, 21), col = c(, , , , , )) abline(lm(log2(Talin.KD.Data\(Mean.Clustered.b5.per.Cell) ~ Talin.KD.Data\)Total.Talin.per.Cell), col = ) r lm.b5.intensity_vs_talin.per.cell <- lm(log2(Talin.KD.Data\(Mean.Clustered.b5.per.Cell) ~ Talin.KD.Data\)Total.Talin.per.Cell) TalinKD_b5_Intensity_panel_file_name <- 3I TalinKD_vs_b5_Intensity
svglite(paste(panel_output_path, TalinKD_b5_Intensity_panel_file_name, .svg, sep = \), width=3, height=3) plot(Talin.KD.Data\(Total.Talin.per.Cell, log2(Talin.KD.Data\)Mean.Clustered.b5.per.Cell), xlim = c(0,1), ylim = c(0,4), pch = c(0, 0, 0, 0, 2, 21), col = c(, , , , , ), cex.axis = 0.5, cex.label = 0.5) r abline(lm(log2(Talin.KD.Data\(Mean.Clustered.b5.per.Cell) ~ Talin.KD.Data\)Total.Talin.per.Cell), col = ) dev.off()
quartz_off_screen
2
rr names(CytoD2h_vs_control) <- c(_Intensity, 5_Intensity, ) Vinculin_Intensity_mu_test <- wilcox.test(Vinculin_Intensity ~ Condition, data = CytoD2h_vs_control)$p.value print(Vinculin_Intensity_mu_test)
[1] 1.785123e-166
rr b5_Intensity_mu_test <- wilcox.test(b5_Intensity ~ Condition, data = CytoD2h_vs_control)$p.value print(b5_Intensity_mu_test)
[1] 1.624755e-99
rr # Set your path to input data data_input_path <- /Users/johnlock/Dropbox/A-CMAC Manchester collaboration/NCB submission documents/Resubmission 2/Source_Data
data_file <- Table 1 Statistics Source Data.xlsx
# Import individual data files Primary.siRNA.Screen.Data <- read_xlsx(paste(data_input_path, data_file, sep = /), sheet = 5B_C, col_names = TRUE) Neo.LY.Drug.Treatment.Data <- read_xlsx(paste(data_input_path, data_file, sep = /), sheet = 5D, col_names = TRUE)
rr siRNA_boxplot_summary_filename <- 5B siRNA boxplot summary
cell.numbers.primary.screen <- count(Primary.siRNA.Screen.Data, siRNA_Target) count(Primary.siRNA.Screen.Data, siRNA_Target) r mean.cell.numbers.primary.screen <- summarise(cell.numbers.primary.screen, avg = mean(n), SD = sd(n), sum(n)) summarise(cell.numbers.primary.screen, avg = mean(n), SD = sd(n), sum(n))
rr ggplot(Primary.siRNA.Screen.Data, aes(x=siRNA_Target, y=Ratio_A.CMAC.to.T.CMAC_mean_intensity.z.score, color = factor(siRNA_Treatment_Type))) + geom_boxplot(notch = TRUE, outlier.colour = ) + scale_x_discrete(limits=c(,,_3Ctrl_Pool, -TARGETplus Non-targeting Control, Non-targeting Control, , 5, , 4KA, 4K2B, , 5K1B, 5K1C, 3C2A)) + theme(axis.text.x = element_text(angle=90)) + ylim(-3, 6)#+ svglite(paste(panel_output_path, siRNA_boxplot_summary_filename, .svg, sep = \), width=7, height=5) r ggplot(Primary.siRNA.Screen.Data, aes(x=siRNA_Target, y=Ratio_A.CMAC.to.T.CMAC_mean_intensity.z.score, color = factor(siRNA_Treatment_Type))) + geom_boxplot(notch = TRUE, outlier.colour = ) + scale_x_discrete(limits=c(,,_3Ctrl_Pool, -TARGETplus Non-targeting Control, Non-targeting Control, , 5, , 4KA, 4K2B, , 5K1B, 5K1C, 3C2A)) + theme(axis.text.x = element_text(angle=90)) + ylim(-3, 6)#+ dev.off()
quartz_off_screen
2
rr siRNA_parcoord_summary_filename <- 5C siRNA parcoord summary
Primary_screen_siRNA_subset <- subset(Primary.siRNA.Screen.Data, siRNA_Target == | siRNA_Target == | siRNA_Target == _3Ctrl_Pool | siRNA_Target == -TARGETplus Non-targeting Control | siRNA_Target == Non-targeting Control | siRNA_Target == | siRNA_Target == 4KA | siRNA_Target == 4K2B | siRNA_Target == | siRNA_Target == 5K1B | siRNA_Target == 5K1C | siRNA_Target == 3C2A) Primary_screen_siRNA_subset_grouped <- aggregate(Primary_screen_siRNA_subset[, c(11,13,18)], list(Primary_screen_siRNA_subset\(siRNA_Target), median) Primary_screen_siRNA_subset_grouped\)Group.1 <- as.factor(Primary_screen_siRNA_subset_grouped$Group.1) names(Primary_screen_siRNA_subset_grouped) <- c(_Target, _Reticular_Adhesion_Intensity, _Focal_Adhesion_Intensity, _Reticular_Focal_Adhesion_Intensity) svglite(paste(panel_output_path, siRNA_parcoord_summary_filename, .svg, sep = \), width=7, height=5) ggparcoord(Primary_screen_siRNA_subset_grouped, columns = c(2,3,4), groupColumn = 1, scale = , order = c(2,3,4)) + theme(axis.text.x = element_text(angle=60, hjust = 1, size = 5)) dev.off()
null device
1
rr ggparcoord(Primary_screen_siRNA_subset_grouped, columns = c(2,3,4), groupColumn = 1, scale = , order = c(2,3,4)) + theme(axis.text.x = element_text(angle=60, hjust = 1, size = 5))
rr siRNA.conditions.for.U.test <- subset(Primary.siRNA.Screen.Data, siRNA_Target != & siRNA_Target != & siRNA_Target != _3Ctrl_Pool & siRNA_Target != Non-targeting Control) p_values <- c() for (i in 1:length(unique(siRNA.conditions.for.U.test\(siRNA_Target))){ for (j in 1:length(unique(siRNA.conditions.for.U.test\)siRNA_Target))){ siRNA_Target_1 = subset(siRNA.conditions.for.U.test, siRNA.conditions.for.U.test\(siRNA_Target == unique(siRNA.conditions.for.U.test\)siRNA_Target)[i]) siRNA_Target_1\(siRNA_Target_Number = as.factor(i) siRNA_Target_2 = subset(siRNA.conditions.for.U.test, siRNA.conditions.for.U.test\)siRNA_Target == unique(siRNA.conditions.for.U.test\(siRNA_Target)[j]) siRNA_Target_2\)siRNA_Target_Number = as.factor(j+100) siRNA_Target_Pair = rbind(siRNA_Target_1, siRNA_Target_2) siRNA_Target_Pair_ttest = wilcox.test(Ratio_A.CMAC.to.T.CMAC_mean_intensity.z.score ~ siRNA_Target_Number, data = siRNA_Target_Pair)\(p.value p_values = rbind(p_values, siRNA_Target_Pair_ttest) } } p_values <- as.vector(p_values) p_value_matrix <- matrix(p_values,nrow = length(unique(siRNA.conditions.for.U.test\)siRNA_Target)),ncol = length(unique(siRNA.conditions.for.U.test\(siRNA_Target))) columns <- as.vector(unique(siRNA.conditions.for.U.test\)siRNA_Target)) rows <- as.vector(unique(siRNA.conditions.for.U.test\(siRNA_Target)) p_value_dataframe <- as.data.frame(p_value_matrix, row.names = rows) names(p_value_dataframe) <- columns # print(p_value_dataframe) # Introduce Holm (also called Holm-Bonferroni) p-value correction (used because it gives same stringency against false positives (type 1 errors) with lower probability of false negatives (type 2 errors))) corrected_p_values <- p.adjust(p_values, method = \holm\, n = length(p_values)) corrected_p_value_matrix2 <- matrix(corrected_p_values,nrow = length(unique(siRNA.conditions.for.U.test\)siRNA_Target)),ncol = length(unique(siRNA.conditions.for.U.test\(siRNA_Target))) columns <- as.vector(unique(siRNA.conditions.for.U.test\)siRNA_Target)) rows <- as.vector(unique(siRNA.conditions.for.U.test$siRNA_Target)) corrected_p_value_dataframe <- as.data.frame(corrected_p_value_matrix2, row.names = rows) names(corrected_p_value_dataframe) <- columns print(corrected_p_value_dataframe)
rr heatmap.2(log(corrected_p_value_matrix2), dendrogram = ‘none’, Rowv = FALSE, Colv = FALSE, labRow = rows, labCol = columns, srtCol = 45, cexRow = 1, cexCol = 1, margins = c(7, 9), symm = TRUE, revC = FALSE, breaks = c(log(1e-300), log(1e-50), log(1e-10), log(1e-6), log(0.001), log(0.01), log(0.05), log(1)), rowsep = 1:nrow(corrected_p_value_matrix2), colsep = 1:ncol(corrected_p_value_matrix2), sepcolor = ‘darkgrey’, sepwidth = c(0.02,0.02), trace = ‘none’, col = c(, , , , , , ), denscol = NULL, keysize = 1.5, key.title = NA)
rr PIP_drug_boxplot_summary_filename <- 5D PIP drug boxplot summary
cell.numbers.NEO_LY.screen <- count(Neo.LY.Drug.Treatment.Data, Drug) count(Neo.LY.Drug.Treatment.Data, Drug) r mean.cell.numbers.NEO_LY.screen <- summarise(cell.numbers.NEO_LY.screen, avg = mean(n), SD = sd(n), sum(n)) summarise(cell.numbers.NEO_LY.screen, avg = mean(n), SD = sd(n), sum(n))
rr svglite(paste(panel_output_path, PIP_drug_boxplot_summary_filename, .svg, sep = \), width=7, height=5) p1 <- ggplot(Neo.LY.Drug.Treatment.Data, aes(x=Drug, y=Median.Mean..z.score, color = factor(Drug))) + geom_boxplot(notch = TRUE, outlier.colour = ) + scale_x_discrete(limits=c(, _10mM,_25mM)) + theme(axis.text.x = element_text(angle=45, vjust = 0.5)) + coord_cartesian(ylim=c(-3, 6)) + ylab(Reticular Adhesion Intensity) + theme(legend.position=) ## plot T-CMAC intensity Z-score ggplot2 boxplot - per Drug p2 <- ggplot(Neo.LY.Drug.Treatment.Data, aes(x=Drug, y=Median.Mean._T.CMACs.z.score, color = factor(Drug))) + geom_boxplot(notch = TRUE, outlier.colour = ) + scale_x_discrete(limits=c(, _10mM,_25mM)) + theme(axis.text.x = element_text(angle=45, vjust = 0.5)) + coord_cartesian(ylim=c(-3, 6)) + ylab(Focal Adhesion Intensity) + theme(legend.position=) ## plot intensity ratio Z-score ggplot2 boxplot - per Drug p3 <- ggplot(Neo.LY.Drug.Treatment.Data, aes(x=Drug, y=Ratio_A.CMAC.to.T.CMAC_mean_intensity.z.score, color = factor(Drug))) + geom_boxplot(notch = TRUE, outlier.colour = ) + scale_x_discrete(limits=c(, _10mM,_25mM)) + theme(axis.text.x = element_text(angle=45, vjust = 0.5)) + coord_cartesian(ylim=c(-3, 6)) + ylab(Reticular:Focal Adhesion Intensity) + theme(legend.position=) grid.arrange(p1, p2, p3, ncol = 3) dev.off()
null device
1
rr p1 <- ggplot(Neo.LY.Drug.Treatment.Data, aes(x=Drug, y=Median.Mean..z.score, color = factor(Drug))) + geom_boxplot(notch = TRUE, outlier.colour = ) + scale_x_discrete(limits=c(, _10mM,_25mM)) + theme(axis.text.x = element_text(angle=45, vjust = 0.5)) + coord_cartesian(ylim=c(-3, 6)) + ylab(Reticular Adhesion Intensity) + theme(legend.position=) ## plot T-CMAC intensity Z-score ggplot2 boxplot - per Drug p2 <- ggplot(Neo.LY.Drug.Treatment.Data, aes(x=Drug, y=Median.Mean._T.CMACs.z.score, color = factor(Drug))) + geom_boxplot(notch = TRUE, outlier.colour = ) + scale_x_discrete(limits=c(, _10mM,_25mM)) + theme(axis.text.x = element_text(angle=45, vjust = 0.5)) + coord_cartesian(ylim=c(-3, 6)) + ylab(Focal Adhesion Intensity) + theme(legend.position=) ## plot intensity ratio Z-score ggplot2 boxplot - per Drug p3 <- ggplot(Neo.LY.Drug.Treatment.Data, aes(x=Drug, y=Ratio_A.CMAC.to.T.CMAC_mean_intensity.z.score, color = factor(Drug))) + geom_boxplot(notch = TRUE, outlier.colour = ) + scale_x_discrete(limits=c(, _10mM,_25mM)) + theme(axis.text.x = element_text(angle=45, vjust = 0.5)) + coord_cartesian(ylim=c(-3, 6)) + ylab(Reticular:Focal Adhesion Intensity) + theme(legend.position=) grid.arrange(p1, p2, p3, ncol = 3)
rr p_values <- c() for (i in 1:length(unique(Neo.LY.Drug.Treatment.Data\(Drug))){ for (j in 1:length(unique(Neo.LY.Drug.Treatment.Data\)Drug))){ Drug_1 = subset(Neo.LY.Drug.Treatment.Data, Neo.LY.Drug.Treatment.Data\(Drug == unique(Neo.LY.Drug.Treatment.Data\)Drug)[i]) Drug_1\(Drug_Number = as.factor(i) Drug_2 = subset(Neo.LY.Drug.Treatment.Data, Neo.LY.Drug.Treatment.Data\)Drug == unique(Neo.LY.Drug.Treatment.Data\(Drug)[j]) Drug_2\)Drug_Number = as.factor(j+100) Drug_Pair = rbind(Drug_1, Drug_2) Drug_Pair_ttest = wilcox.test(Median.Mean..z.score ~ Drug_Number, data = Drug_Pair)\(p.value p_values = rbind(p_values, Drug_Pair_ttest) } } p_values <- as.vector(p_values) p_value_matrix <- matrix(p_values,nrow = length(unique(Neo.LY.Drug.Treatment.Data\)Drug)),ncol = length(unique(Neo.LY.Drug.Treatment.Data\(Drug))) columns <- as.vector(unique(Neo.LY.Drug.Treatment.Data\)Drug)) rows <- as.vector(unique(Neo.LY.Drug.Treatment.Data\(Drug)) p_value_dataframe <- as.data.frame(p_value_matrix, row.names = rows) names(p_value_dataframe) <- columns # print(p_value_dataframe) # Introduce Holm (also called Holm-Bonferroni) p-value correction (used because it gives same stringency against false positives (type 1 errors) with lower probability of false negatives (type 2 errors))) corrected_p_values <- p.adjust(p_values, method = \holm\, n = length(p_values)) corrected_p_value_matrix3 <- matrix(corrected_p_values,nrow = length(unique(Neo.LY.Drug.Treatment.Data\)Drug)),ncol = length(unique(Neo.LY.Drug.Treatment.Data\(Drug))) columns <- as.vector(unique(Neo.LY.Drug.Treatment.Data\)Drug)) rows <- as.vector(unique(Neo.LY.Drug.Treatment.Data$Drug)) corrected_p_value_dataframe <- as.data.frame(corrected_p_value_matrix3, row.names = rows) names(corrected_p_value_dataframe) <- columns print(corrected_p_value_dataframe)
rr p_values <- c() for (i in 1:length(unique(Neo.LY.Drug.Treatment.Data\(Drug))){ for (j in 1:length(unique(Neo.LY.Drug.Treatment.Data\)Drug))){ Drug_1 = subset(Neo.LY.Drug.Treatment.Data, Neo.LY.Drug.Treatment.Data\(Drug == unique(Neo.LY.Drug.Treatment.Data\)Drug)[i]) Drug_1\(Drug_Number = as.factor(i) Drug_2 = subset(Neo.LY.Drug.Treatment.Data, Neo.LY.Drug.Treatment.Data\)Drug == unique(Neo.LY.Drug.Treatment.Data\(Drug)[j]) Drug_2\)Drug_Number = as.factor(j+100) Drug_Pair = rbind(Drug_1, Drug_2) Drug_Pair_ttest = wilcox.test(Median.Mean._T.CMACs.z.score ~ Drug_Number, data = Drug_Pair)\(p.value p_values = rbind(p_values, Drug_Pair_ttest) } } p_values <- as.vector(p_values) p_value_matrix <- matrix(p_values,nrow = length(unique(Neo.LY.Drug.Treatment.Data\)Drug)),ncol = length(unique(Neo.LY.Drug.Treatment.Data\(Drug))) columns <- as.vector(unique(Neo.LY.Drug.Treatment.Data\)Drug)) rows <- as.vector(unique(Neo.LY.Drug.Treatment.Data\(Drug)) p_value_dataframe <- as.data.frame(p_value_matrix, row.names = rows) names(p_value_dataframe) <- columns # print(p_value_dataframe) # Introduce Holm (also called Holm-Bonferroni) p-value correction (used because it gives same stringency against false positives (type 1 errors) with lower probability of false negatives (type 2 errors))) corrected_p_values <- p.adjust(p_values, method = \holm\, n = length(p_values)) corrected_p_value_matrix4 <- matrix(corrected_p_values,nrow = length(unique(Neo.LY.Drug.Treatment.Data\)Drug)),ncol = length(unique(Neo.LY.Drug.Treatment.Data\(Drug))) columns <- as.vector(unique(Neo.LY.Drug.Treatment.Data\)Drug)) rows <- as.vector(unique(Neo.LY.Drug.Treatment.Data$Drug)) corrected_p_value_dataframe <- as.data.frame(corrected_p_value_matrix4, row.names = rows) names(corrected_p_value_dataframe) <- columns print(corrected_p_value_dataframe)
rr p_values <- c() for (i in 1:length(unique(Neo.LY.Drug.Treatment.Data\(Drug))){ for (j in 1:length(unique(Neo.LY.Drug.Treatment.Data\)Drug))){ Drug_1 = subset(Neo.LY.Drug.Treatment.Data, Neo.LY.Drug.Treatment.Data\(Drug == unique(Neo.LY.Drug.Treatment.Data\)Drug)[i]) Drug_1\(Drug_Number = as.factor(i) Drug_2 = subset(Neo.LY.Drug.Treatment.Data, Neo.LY.Drug.Treatment.Data\)Drug == unique(Neo.LY.Drug.Treatment.Data\(Drug)[j]) Drug_2\)Drug_Number = as.factor(j+100) Drug_Pair = rbind(Drug_1, Drug_2) Drug_Pair_ttest = wilcox.test(as.numeric(Ratio_A.CMAC.to.T.CMAC_mean_intensity.z.score) ~ Drug_Number, data = Drug_Pair)\(p.value p_values = rbind(p_values, Drug_Pair_ttest) } } p_values <- as.vector(p_values) p_value_matrix <- matrix(p_values,nrow = length(unique(Neo.LY.Drug.Treatment.Data\)Drug)),ncol = length(unique(Neo.LY.Drug.Treatment.Data\(Drug))) columns <- as.vector(unique(Neo.LY.Drug.Treatment.Data\)Drug)) rows <- as.vector(unique(Neo.LY.Drug.Treatment.Data\(Drug)) p_value_dataframe <- as.data.frame(p_value_matrix, row.names = rows) names(p_value_dataframe) <- columns # print(p_value_dataframe) # Introduce Holm (also called Holm-Bonferroni) p-value correction (used because it gives same stringency against false positives (type 1 errors) with lower probability of false negatives (type 2 errors))) corrected_p_values <- p.adjust(p_values, method = \holm\, n = length(p_values)) corrected_p_value_matrix5 <- matrix(corrected_p_values,nrow = length(unique(Neo.LY.Drug.Treatment.Data\)Drug)),ncol = length(unique(Neo.LY.Drug.Treatment.Data\(Drug))) columns <- as.vector(unique(Neo.LY.Drug.Treatment.Data\)Drug)) rows <- as.vector(unique(Neo.LY.Drug.Treatment.Data$Drug)) corrected_p_value_dataframe <- as.data.frame(corrected_p_value_matrix5, row.names = rows) names(corrected_p_value_dataframe) <- columns print(corrected_p_value_dataframe)
rr # Set your path to input data data_input_path <- /Users/johnlock/Dropbox/A-CMAC Manchester collaboration/NCB submission documents/Resubmission 2/Source_Data
data_file <- Table 1 Statistics Source Data.xlsx
# Import individual data files Proliferation <- read_xlsx(paste(data_input_path, data_file, sep = /), sheet = 6A, col_names = TRUE) Proliferation\(Day <- as.factor(Proliferation\)Day) EDU_Incorporation_Percentage <- read_xlsx(paste(data_input_path, data_file, sep = /), sheet = 6B, col_names = TRUE) EDU_Incorporation_Percentage <- EDU_Incorporation_Percentage[,1:2]
rr Prol_Sum <- Proliferation %>% group_by(.dots = c(‘Condition’, ‘Day’)) %>% summarise( mean_proliferation = mean(Relative.Proliferation, na.rm = TRUE), SD_proliferation = sd(Relative.Proliferation, na.rm = TRUE) )
rr proliferation_linechart_filename <- 6A proliferation linechart
# The errorbars overlapped, so use position_dodge to move them horizontally pd <- position_dodge(0.1) # move them .05 to the left and right svglite(paste(panel_output_path, proliferation_linechart_filename, .svg, sep = \), width=4, height=3) ggplot(Prol_Sum, aes(x=Day, y=mean_proliferation, colour=Condition, group=Condition)) + geom_errorbar(aes(ymin=mean_proliferation-SD_proliferation, ymax=mean_proliferation+SD_proliferation), width=.1, position=pd) + scale_color_manual(values = c( = #49494A, 5KD = #EC1C24)) + geom_line(position=pd, size = 1) + geom_point(position=pd, size=3, shape=19) + # 21 is filled circle xlab() + ylab(Proliferation) + theme_light() dev.off()
null device
1
rr ggplot(Prol_Sum, aes(x=Day, y=mean_proliferation, colour=Condition, group=Condition)) + geom_errorbar(aes(ymin=mean_proliferation-SD_proliferation, ymax=mean_proliferation+SD_proliferation), width=.1, position=pd) + scale_color_manual(values = c( = #49494A, 5KD = #EC1C24)) + geom_line(position=pd, size = 1) + geom_point(position=pd, size=3, shape=19) + # 21 is filled circle xlab() + ylab(Proliferation) + theme_light()
rr Proliferation\(Day <- as.factor(Proliferation\)Day) p_values <- c() for (i in 1:length(unique(Proliferation$Day))){
single_day_data = subset(Proliferation, Proliferation$Day == unique(Proliferation$Day)[i])
Control_vs_b5KD_per_day_ttest = t.test(Relative.Proliferation ~ Condition, data = single_day_data, paired = FALSE)$p.value
p_values = rbind(p_values, Control_vs_b5KD_per_day_ttest)
} p_values <- as.vector(p_values) p_value_matrix <- matrix(p_values,nrow = 1,ncol = length(p_values)) columns <- as.vector(paste(, unique(Proliferation\(Day), sep = \ \)) rows <- paste(unique(Proliferation\)Condition)[1], unique(Proliferation$Condition)[2], sep = _vs_) p_value_dataframe <- as.data.frame(p_value_matrix, row.names = rows) names(p_value_dataframe) <- columns print(p_value_dataframe)
rr #define panel_file_name EDU_boxplot_panel_file_name <- 6B EDU Incorporation Boxplots
EDU_Incorporation_Percentage\(Condition <- factor(EDU_Incorporation_Percentage\)Condition, c(, 5KD)) svglite(paste(panel_output_path, EDU_boxplot_panel_file_name, .svg, sep = \), width=4, height=6) boxplot(Percentage.EDU.positive.Cells.per.image ~ Condition, data = EDU_Incorporation_Percentage, notch = T, varwidth = T, outline = F, ylab = % EDU positive cells, col = c(control = #49494A, b5KD = #EC1C24), ylim = c(0,115)) beeswarm(Percentage.EDU.positive.Cells.per.image ~ Condition, data = EDU_Incorporation_Percentage, method = ‘center’, add = T, col = 2, pch = 1, cex = 1.3) dev.off()
null device
1
rr boxplot(Percentage.EDU.positive.Cells.per.image ~ Condition, data = EDU_Incorporation_Percentage, notch = T, varwidth = T, outline = F, ylab = % EDU positive cells, col = c(control = #49494A, b5KD = #EC1C24), ylim = c(0,115)) r beeswarm(Percentage.EDU.positive.Cells.per.image ~ Condition, data = EDU_Incorporation_Percentage, method = ‘center’, add = T, col = 2, pch = 1, cex = 1.3)
rr #Note: Equivalent to Wilcoxon rank sum test as data is unpaired Utest.b5KD.vs.control_EDU <- wilcox.test(Percentage.EDU.positive.Cells.per.image ~ Condition, data = EDU_Incorporation_Percentage)
cannot compute exact p-value with ties
rr Mann_Whitney_EDU_Incorporation_Score <- Utest.b5KD.vs.control_EDU$p.value print(Mann_Whitney_EDU_Incorporation_Score)
[1] 0.5337576
rr # Set your path to input data data_input_path <- /Users/johnlock/Dropbox/A-CMAC Manchester collaboration/NCB submission documents/Resubmission 2/Source_Data
data_file <- Table 1 Statistics Source Data.xlsx
# Import individual data files STORM.Data <- read_xlsx(paste(data_input_path, data_file, sep = /), sheet = 7E_F, col_names = TRUE)
rr ## Determine CMAC.type (Reticular vs Focal depending on presence of A in \(File (column)) CMAC.type = c() for (cluster in 1:length(STORM.Data\)File) ) { CMAC.type[cluster] <- if(grepl(, STORM.Data\(File[cluster]) == TRUE) { \Reticular\
} else if(grepl(\T\, STORM.Data\)File[cluster]) == TRUE) {
} else if(grepl(, STORM.Data\(File[cluster]) == TRUE) { \Non-Retraction\
} else { \Retraction\
} } STORM.Data\)CMAC.type <- cbind(as.character(CMAC.type)) ## Unique nanocluster ID generation = merge of Folder, File and Cluster STORM.Data\(Unique.ncID <- as.factor(paste(STORM.Data\)Folder, STORM.Data\(File, STORM.Data\)Cluster)) STORM.Data\(Unique.CMACID <- as.factor(paste(STORM.Data\)Folder, STORM.Data$File))
rr ## NND analysis of Incite3 data NNND <- c() for (i in 1:length(unique(STORM.Data\(Unique.CMACID)) ) { temp = STORM.Data[STORM.Data\)Unique.CMACID == unique(STORM.Data\(Unique.CMACID)[i], ] NNND = append(NNND, nndist(temp\)Centroid.nm., temp$Centroid.nm..1, k=1)) }
NAs introduced by coercionNAs introduced by coercion
rr STORM.Data$NNND <- cbind(NNND)
rr NNND.A = subset(STORM.Data\(NNND, CMAC.type == \Non-Retraction\) NNND.T = subset(STORM.Data\)NNND, CMAC.type == ) TEST_nnnd_NRvsR <- wilcox.test(NNND.A, NNND.T, paired = FALSE) NNND_boxplot_filename <- 7E NNND boxplot
svglite(paste(panel_output_path, NNND_boxplot_filename, .svg, sep = \), width=4, height=3) boxplot(NNND ~ CMAC.type, subset(STORM.Data, CMAC.type == -Retraction | CMAC.type == ), notch = TRUE, outline = FALSE, main = Neighbour Distance of Nanocluster by CMAC Type, sub = TEST_nnnd_NRvsR$p.value, xlab = type, ylab = Neighbour Distance (nm), plot =
) dev.off()
null device
1
rr boxplot(NNND ~ CMAC.type, subset(STORM.Data, CMAC.type == -Retraction | CMAC.type == ), notch = TRUE, outline = FALSE, main = Neighbour Distance of Nanocluster by CMAC Type, sub = TEST_nnnd_NRvsR$p.value, xlab = type, ylab = Neighbour Distance (nm), plot =
)
rr STORM.Data\(Molecules <- as.numeric(STORM.Data\)Molecules)
NAs introduced by coercion
rr Molecules.A = subset(STORM.Data\(Molecules, CMAC.type == \Non-Retraction\) Molecules.T = subset(STORM.Data\)Molecules, CMAC.type == ) Mol_Count_boxplot_filename <- 7F Molecule Count boxplot
svglite(paste(panel_output_path, Mol_Count_boxplot_filename, .svg, sep = \), width=4, height=3) TEST_Molecules.A_RetvsFoc <- wilcox.test(Molecules.A, Molecules.T, paired = FALSE) boxplot(Molecules ~ CMAC.type, subset(STORM.Data, CMAC.type == -Retraction | CMAC.type == ), notch = TRUE,outline = FALSE, main = of Integrin avb5 Molecules by CMAC Type, sub = TEST_Molecules.A_RetvsFoc$p.value, xlab = type, ylab = of Integrin avb5 Molecules, plot =
) dev.off()
null device
1
rr TEST_Molecules.A_RetvsFoc <- wilcox.test(Molecules.A, Molecules.T, paired = FALSE) boxplot(Molecules ~ CMAC.type, subset(STORM.Data, CMAC.type == -Retraction | CMAC.type == ), notch = TRUE,outline = FALSE, main = of Integrin avb5 Molecules by CMAC Type, sub = TEST_Molecules.A_RetvsFoc$p.value, xlab = type, ylab = of Integrin avb5 Molecules, plot =
)
rr # Set your path to input data data_input_path <- /Users/johnlock/Dropbox/A-CMAC Manchester collaboration/NCB submission documents/Resubmission 2/Source_Data
data_file <- Table 1 Statistics Source Data.xlsx
# Import individual data files Residual_Angle_data <- read_xlsx(paste(data_input_path, data_file, sep = /), sheet = 8A_B, col_names = TRUE) Cell_Division_Defect_Quantification <- read_xlsx(paste(data_input_path, data_file, sep = /), sheet = 8C, col_names = TRUE) Cell_Division_Defect_Quantification <- Cell_Division_Defect_Quantification[1:9, 1:4]
rr #define panel_file_name boxplot_panel_file_name <- 8A Residual Angle Boxplots
# oldpar <- par() # par(mfrow = c(1,1), mar = c(5,5,1,1)) Residual_Angle_dataCN <- Residual_Angle_data # Set Condition order Residual_Angle_dataCN\(Condition <- factor(Residual_Angle_dataCN\)Condition, c(, 5KD, )) svglite(paste(panel_output_path, boxplot_panel_file_name, .svg, sep = \), width=4, height=6) boxplot(Pre.to.Post.Mitosis.Angle ~ Condition, data = Residual_Angle_dataCN, notch = T, varwidth = T, outline = F, ylab = Angle, col = c(control = #49494A, b5KD = #EC1C24, Rescue = #BCBCBC), ylim = c(0,115)) beeswarm(Pre.to.Post.Mitosis.Angle ~ Condition, data = Residual_Angle_dataCN, method = ‘center’, add = T, col = 2, pch = 1, cex = 0.5) dev.off()
null device
1
rr boxplot(Pre.to.Post.Mitosis.Angle ~ Condition, data = Residual_Angle_dataCN, notch = T, varwidth = T, outline = F, ylab = Angle, col = c(control = #49494A, b5KD = #EC1C24, Rescue = #BCBCBC), ylim = c(0,115)) r beeswarm(Pre.to.Post.Mitosis.Angle ~ Condition, data = Residual_Angle_dataCN, method = ‘center’, add = T, col = 2, pch = 1, cex = 0.5)
rr #define panel_file_name density_panel_file_name <- 8B Residual Angle Density Plot
svglite(paste(panel_output_path, density_panel_file_name, .svg, sep = \), width=4, height=6) sm.density.compare(Residual_Angle_dataCN\(Pre.to.Post.Mitosis.Angle, Residual_Angle_dataCN\)Condition, col = c(#49494A, #EC1C24, #BCBCBC), lty = c(1,1,1), lwd = 2, h = 5, xlim = c(0,90), xlab = Angle) colfill<-c(#49494A, #EC1C24, #BCBCBC) legend(x=40, y=0.033, levels(Residual_Angle_dataCN$Condition), fill=colfill, bty = ) dev.off()
null device
1
rr sm.density.compare(Residual_Angle_dataCN\(Pre.to.Post.Mitosis.Angle, Residual_Angle_dataCN\)Condition, col = c(#49494A, #EC1C24, #BCBCBC), lty = c(1,1,1), lwd = 2, h = 5, xlim = c(0,90), xlab = Angle) r colfill<-c(#49494A, #EC1C24, #BCBCBC) legend(x=40, y=0.033, levels(Residual_Angle_dataCN$Condition), fill=colfill, bty = )
rr Residual_Angle_dataCN %>% group_by(Condition) %>% summarise( n = n())
rr #Note: Equivalent to Wilcoxon rank sum test as data is unpaired b5KD_Control <- subset(Residual_Angle_data, Residual_Angle_data\(Condition == \b5KD\ | Residual_Angle_data\)Condition == ) Utest.b5KD.vs.control <- wilcox.test(Pre.to.Post.Mitosis.Angle ~ Condition, data = b5KD_Control) Utest.b5KD.vs.control.pval <- Utest.b5KD.vs.control\(p.value b5KD_Rescue <- subset(Residual_Angle_data, Residual_Angle_data\)Condition == 5KD | Residual_Angle_data\(Condition == \Rescue\) Utest.b5KD.vs.Rescue <- wilcox.test(Pre.to.Post.Mitosis.Angle ~ Condition, data = b5KD_Rescue) Utest.b5KD.vs.Rescue.pval <- Utest.b5KD.vs.Rescue\)p.value Control_Rescue <- subset(Residual_Angle_data, Residual_Angle_data\(Condition == \Control\ | Residual_Angle_data\)Condition == ) Utest.control.vs.Rescue <- wilcox.test(Pre.to.Post.Mitosis.Angle ~ Condition, data = Control_Rescue) Utest.control.vs.Rescue.pval <- Utest.control.vs.Rescue\(p.value Mann_Whitney_Scores <- data.frame(Condition_Comparison = c(\b5KD_vs_Control\, \b5KD_vs_Rescue\, \Control_vs_Rescue\)) Mann_Whitney_Scores\)pvalue <- rbind(Utest.b5KD.vs.control.pval, Utest.b5KD.vs.Rescue.pval, Utest.control.vs.Rescue.pval) print(Mann_Whitney_Scores)
rr Residual_Angle_dataFN <- Residual_Angle_data Residual_Angle_dataFN\(File.Name <- factor(Residual_Angle_dataFN\)File.Name, c(_3_MMStack_control Hela 01.ome.tif, _3_MMStack_control Hela 02.ome.tif, _3_MMStack_b5 KD plus GFP 01.ome.tif, _3_MMStack_b5 KD plus GFP 02.ome.tif, _3_MMStack_b5 KD plus WTb5-GFP 01.ome.tif, _3_MMStack_b5 KD plus WTb5-GFP 02.ome.tif)) boxplot(Pre.to.Post.Mitosis.Angle ~ File.Name, data = Residual_Angle_dataFN, notch = T, varwidth = T, outline = F, ylab = Angle, col = c(#49494A, #49494A, #EC1C24, #EC1C24,#BCBCBC, #BCBCBC), ylim = c(0,115), names = c(_01, _02, 5KD_01, 5KD_02, _01, _02))
rr # beeswarm(Pre.to.Post.Mitosis.Angle ~ File.Name, data = Residual_Angle_dataFN, method = ‘swarm’, add = T, col = 2, pch = 0, cex = 0.7)
rr sm.density.compare(Residual_Angle_dataFN\(Pre.to.Post.Mitosis.Angle, Residual_Angle_dataFN\)File.Name, col = c(#49494A, #49494A, #EC1C24, #EC1C24,#BCBCBC,#BCBCBC), h=5, lty = c(1,5,1,5,1,5), lwd = 2, xlim = c(0,90), xlab = Angle) legend(x=23, y= 0.05, levels(Residual_Angle_dataFN$File.Name), fill=c(#49494A, #49494A, #EC1C24, #EC1C24,#BCBCBC,#BCBCBC), bty = ‘n’)
rr Residual_Angle_dataFN %>% group_by(File.Name) %>% summarise( n = n())
rr CDD_boxplot_panel_file_name <- 8C Cell Division Defect Rate Boxplot
Cell_Division_Defect_Quantification\(Condition <- factor(Cell_Division_Defect_Quantification\)Condition, c(, 5KD, )) svglite(paste(panel_output_path, CDD_boxplot_panel_file_name, .svg, sep = \), width=4, height=6) boxplot(normal.cell.division ~ Condition, data = Cell_Division_Defect_Quantification, outline = FALSE, boxlty = 0, whisklty = 0, staplelty = 0, ylab = Division Percentage, ylim = c(0,115), medlwd = 0.0) beeswarm(normal.cell.division ~ Condition, data = Cell_Division_Defect_Quantification, method = ‘center’, add = T, col = c(control = #49494A, b5KD = #EC1C24, Rescue = #BCBCBC), ylim = c(0,115), pch = 19, cex = 1.3) dev.off()
null device
1
rr boxplot(normal.cell.division ~ Condition, data = Cell_Division_Defect_Quantification, outline = FALSE, boxlty = 0, whisklty = 0, staplelty = 0, ylab = Division Percentage, ylim = c(0,115), medlwd = 0.0) r beeswarm(normal.cell.division ~ Condition, data = Cell_Division_Defect_Quantification, method = ‘center’, add = T, col = c(control = #49494A, b5KD = #EC1C24, Rescue = #BCBCBC), ylim = c(0,115), pch = 19, cex = 1.3)
rr Division_Defect_Summary <- Cell_Division_Defect_Quantification %>% group_by(Condition) %>% summarise( Normal_Division = mean(normal.cell.division), Abnormal_Division = mean(abnormal.cell.division) )
print(Division_Defect_Summary)
rr b5KD_Control_CDD <- subset(Cell_Division_Defect_Quantification, Condition == 5KD | Condition == ) T_test_b5KD_Control_CDD <- t.test(normal.cell.division ~ Condition, data = b5KD_Control_CDD)\(p.value b5KD_Rescue_CDD <- subset(Cell_Division_Defect_Quantification, Condition == \b5KD\ | Condition == \Rescue\) T_test_b5KD_Rescue_CDD <- t.test(normal.cell.division ~ Condition, data = b5KD_Rescue_CDD)\)p.value Control_Rescue_CDD <- subset(Cell_Division_Defect_Quantification, Condition == | Condition == ) T_test_Control_Rescue_CDD <- t.test(normal.cell.division ~ Condition, data = Control_Rescue_CDD)\(p.value Ttest_Scores_CDD <- data.frame(Condition_Comparison = c(\b5KD_Control_CDD\, \b5KD_Rescue_CDD\, \Control_Rescue_CDD\)) Ttest_Scores_CDD\)pvalue <- rbind(T_test_b5KD_Control_CDD, T_test_b5KD_Rescue_CDD, T_test_Control_Rescue_CDD) print(Ttest_Scores_CDD)
rr # Set your path to input data data_input_path <- /Users/johnlock/Dropbox/A-CMAC Manchester collaboration/NCB submission documents/Resubmission 2/Source_Data
data_file <- Table 1 Statistics Source Data.xlsx
# Import individual data files Arp23_Inhibition_Data <- read_xlsx(paste(data_input_path, data_file, sep = /), sheet = 5C, col_names = TRUE) Arp23_Inhibition_Data\(Drug <- recode(Arp23_Inhibition_Data\)Drug, Arp_Inhib_CTRL_689 = _CTRL_689)
rr Arp23_b5_intensity_boxplots_filename <- SF5C Arp23_b5_intensity_boxplots
cell.numbers.ARP23.screen <- count(Arp23_Inhibition_Data, Drug) count(Arp23_Inhibition_Data, Drug) r mean.cell.numbers.ARP23.screen <- summarise(cell.numbers.ARP23.screen, avg = mean(n), SD = sd(n), sum(n)) summarise(cell.numbers.ARP23.screen, avg = mean(n), SD = sd(n), sum(n))
rr svglite(paste(panel_output_path, Arp23_b5_intensity_boxplots_filename, .svg, sep = \), width=7, height=5) p1 <- ggplot(Arp23_Inhibition_Data, aes(x=Drug, y=Median.Mean..z.score, color = factor(Drug))) + geom_boxplot(notch = TRUE, outlier.colour = ) + scale_x_discrete(limits=c(_CTRL_689, _Inhib_666)) + theme(axis.text.x = element_text(angle=45, hjust = 1)) + ylim(-3, 6) + ylab(Reticular Adhesion Intensity) + theme(legend.position=) ## plot T-CMAC intensity Z-score ggplot2 boxplot - per Drug p2 <- ggplot(Arp23_Inhibition_Data, aes(x=Drug, y=Median.Mean._T.CMACs.z.score, color = factor(Drug))) + geom_boxplot(notch = TRUE, outlier.colour = ) + scale_x_discrete(limits=c(_CTRL_689, _Inhib_666)) + theme(axis.text.x = element_text(angle=45, hjust = 1)) + ylim(-3, 6) + ylab(Focal Adhesion Intensity) + theme(legend.position=) ## plot intensity ratio Z-score ggplot2 boxplot - per Drug p3 <- ggplot(Arp23_Inhibition_Data, aes(x=Drug, y=Ratio_A.CMAC.to.T.CMAC_mean_intensity.z.score, color = factor(Drug))) + geom_boxplot(notch = TRUE, outlier.colour = ) + scale_x_discrete(limits=c(_CTRL_689, _Inhib_666)) + theme(axis.text.x = element_text(angle=45, hjust = 1)) + ylim(-3, 6) + ylab(Reticular:Focal Adhesion Intensity) + theme(legend.position=) grid.arrange(p1, p2, p3, ncol = 3) dev.off()
null device
1
rr p1 <- ggplot(Arp23_Inhibition_Data, aes(x=Drug, y=Median.Mean..z.score, color = factor(Drug))) + geom_boxplot(notch = TRUE, outlier.colour = ) + scale_x_discrete(limits=c(_CTRL_689, _Inhib_666)) + theme(axis.text.x = element_text(angle=45, hjust = 1)) + ylim(-3, 6) + ylab(Reticular Adhesion Intensity) + theme(legend.position=) ## plot T-CMAC intensity Z-score ggplot2 boxplot - per Drug p2 <- ggplot(Arp23_Inhibition_Data, aes(x=Drug, y=Median.Mean._T.CMACs.z.score, color = factor(Drug))) + geom_boxplot(notch = TRUE, outlier.colour = ) + scale_x_discrete(limits=c(_CTRL_689, _Inhib_666)) + theme(axis.text.x = element_text(angle=45, hjust = 1)) + ylim(-3, 6) + ylab(Focal Adhesion Intensity) + theme(legend.position=) ## plot intensity ratio Z-score ggplot2 boxplot - per Drug p3 <- ggplot(Arp23_Inhibition_Data, aes(x=Drug, y=Ratio_A.CMAC.to.T.CMAC_mean_intensity.z.score, color = factor(Drug))) + geom_boxplot(notch = TRUE, outlier.colour = ) + scale_x_discrete(limits=c(_CTRL_689, _Inhib_666)) + theme(axis.text.x = element_text(angle=45, hjust = 1)) + ylim(-3, 6) + ylab(Reticular:Focal Adhesion Intensity) + theme(legend.position=) grid.arrange(p1, p2, p3, ncol = 3)
rr p_values <- c() for (i in 1:length(unique(Arp23_Inhibition_Data\(Drug))){ for (j in 1:length(unique(Arp23_Inhibition_Data\)Drug))){ Drug_1 = subset(Arp23_Inhibition_Data, Arp23_Inhibition_Data\(Drug == unique(Arp23_Inhibition_Data\)Drug)[i]) Drug_1\(Drug_Number = as.factor(i) Drug_2 = subset(Arp23_Inhibition_Data, Arp23_Inhibition_Data\)Drug == unique(Arp23_Inhibition_Data\(Drug)[j]) Drug_2\)Drug_Number = as.factor(j+100) Drug_Pair = rbind(Drug_1, Drug_2) Drug_Pair_ttest = wilcox.test(Median.Mean..z.score ~ Drug_Number, data = Drug_Pair)\(p.value p_values = rbind(p_values, Drug_Pair_ttest) } } p_values <- as.vector(p_values) p_value_matrix <- matrix(p_values,nrow = length(unique(Arp23_Inhibition_Data\)Drug)),ncol = length(unique(Arp23_Inhibition_Data\(Drug))) columns <- as.vector(unique(Arp23_Inhibition_Data\)Drug)) rows <- as.vector(unique(Arp23_Inhibition_Data\(Drug)) p_value_dataframe <- as.data.frame(p_value_matrix, row.names = rows) names(p_value_dataframe) <- columns # print(p_value_dataframe) # Introduce Holm (also called Holm-Bonferroni) p-value correction (used because it gives same stringency against false positives (type 1 errors) with lower probability of false negatives (type 2 errors))) corrected_p_values <- p.adjust(p_values, method = \holm\, n = length(p_values)) corrected_p_value_matrix6 <- matrix(corrected_p_values,nrow = length(unique(Arp23_Inhibition_Data\)Drug)),ncol = length(unique(Arp23_Inhibition_Data\(Drug))) columns <- as.vector(unique(Arp23_Inhibition_Data\)Drug)) rows <- as.vector(unique(Arp23_Inhibition_Data$Drug)) corrected_p_value_dataframe <- as.data.frame(corrected_p_value_matrix6, row.names = rows) names(corrected_p_value_dataframe) <- columns print(corrected_p_value_dataframe)
rr p_values <- c() for (i in 1:length(unique(Arp23_Inhibition_Data\(Drug))){ for (j in 1:length(unique(Arp23_Inhibition_Data\)Drug))){ Drug_1 = subset(Arp23_Inhibition_Data, Arp23_Inhibition_Data\(Drug == unique(Arp23_Inhibition_Data\)Drug)[i]) Drug_1\(Drug_Number = as.factor(i) Drug_2 = subset(Arp23_Inhibition_Data, Arp23_Inhibition_Data\)Drug == unique(Arp23_Inhibition_Data\(Drug)[j]) Drug_2\)Drug_Number = as.factor(j+100) Drug_Pair = rbind(Drug_1, Drug_2) Drug_Pair_ttest = wilcox.test(Median.Mean._T.CMACs.z.score ~ Drug_Number, data = Drug_Pair)\(p.value p_values = rbind(p_values, Drug_Pair_ttest) } } p_values <- as.vector(p_values) p_value_matrix <- matrix(p_values,nrow = length(unique(Arp23_Inhibition_Data\)Drug)),ncol = length(unique(Arp23_Inhibition_Data\(Drug))) columns <- as.vector(unique(Arp23_Inhibition_Data\)Drug)) rows <- as.vector(unique(Arp23_Inhibition_Data\(Drug)) p_value_dataframe <- as.data.frame(p_value_matrix, row.names = rows) names(p_value_dataframe) <- columns # print(p_value_dataframe) # Introduce Holm (also called Holm-Bonferroni) p-value correction (used because it gives same stringency against false positives (type 1 errors) with lower probability of false negatives (type 2 errors))) corrected_p_values <- p.adjust(p_values, method = \holm\, n = length(p_values)) corrected_p_value_matrix7 <- matrix(corrected_p_values,nrow = length(unique(Arp23_Inhibition_Data\)Drug)),ncol = length(unique(Arp23_Inhibition_Data\(Drug))) columns <- as.vector(unique(Arp23_Inhibition_Data\)Drug)) rows <- as.vector(unique(Arp23_Inhibition_Data$Drug)) corrected_p_value_dataframe <- as.data.frame(corrected_p_value_matrix7, row.names = rows) names(corrected_p_value_dataframe) <- columns print(corrected_p_value_dataframe)
rr p_values <- c() for (i in 1:length(unique(Arp23_Inhibition_Data\(Drug))){ for (j in 1:length(unique(Arp23_Inhibition_Data\)Drug))){ Drug_1 = subset(Arp23_Inhibition_Data, Arp23_Inhibition_Data\(Drug == unique(Arp23_Inhibition_Data\)Drug)[i]) Drug_1\(Drug_Number = as.factor(i) Drug_2 = subset(Arp23_Inhibition_Data, Arp23_Inhibition_Data\)Drug == unique(Arp23_Inhibition_Data\(Drug)[j]) Drug_2\)Drug_Number = as.factor(j+100) Drug_Pair = rbind(Drug_1, Drug_2) Drug_Pair_ttest = wilcox.test(Ratio_A.CMAC.to.T.CMAC_mean_intensity.z.score ~ Drug_Number, data = Drug_Pair)\(p.value p_values = rbind(p_values, Drug_Pair_ttest) } } p_values <- as.vector(p_values) p_value_matrix <- matrix(p_values,nrow = length(unique(Arp23_Inhibition_Data\)Drug)),ncol = length(unique(Arp23_Inhibition_Data\(Drug))) columns <- as.vector(unique(Arp23_Inhibition_Data\)Drug)) rows <- as.vector(unique(Arp23_Inhibition_Data\(Drug)) p_value_dataframe <- as.data.frame(p_value_matrix, row.names = rows) names(p_value_dataframe) <- columns # print(p_value_dataframe) # Introduce Holm (also called Holm-Bonferroni) p-value correction (used because it gives same stringency against false positives (type 1 errors) with lower probability of false negatives (type 2 errors))) corrected_p_values <- p.adjust(p_values, method = \holm\, n = length(p_values)) corrected_p_value_matrix8 <- matrix(corrected_p_values,nrow = length(unique(Arp23_Inhibition_Data\)Drug)),ncol = length(unique(Arp23_Inhibition_Data\(Drug))) columns <- as.vector(unique(Arp23_Inhibition_Data\)Drug)) rows <- as.vector(unique(Arp23_Inhibition_Data$Drug)) corrected_p_value_dataframe <- as.data.frame(corrected_p_value_matrix8, row.names = rows) names(corrected_p_value_dataframe) <- columns print(corrected_p_value_dataframe)
rr # Set your path to input data data_input_path <- /Users/johnlock/Dropbox/A-CMAC Manchester collaboration/NCB submission documents/Resubmission 2/Source_Data
data_file <- Table 1 Statistics Source Data.xlsx
# Import individual data files Secondary.siRNA.Screen.Data <- read_xlsx(paste(data_input_path, data_file, sep = /), sheet = 6A, col_names = TRUE)
rr Secondary_siRNA_Screen_boxplots_filename <- SF6A Secondary_siRNA_Screen_boxplot_summary
Secondary.siRNA.Screen.Data\(siRNA_Target <- recode(Secondary.siRNA.Screen.Data\)siRNA_Target, PIK3C2a_1 = 3C2A_1) Secondary.siRNA.Screen.Data\(siRNA_Target <- recode(Secondary.siRNA.Screen.Data\)siRNA_Target, PIK3C2a_2 = 3C2A_2) Secondary.siRNA.Screen.Data\(siRNA_Target <- recode(Secondary.siRNA.Screen.Data\)siRNA_Target, PIK3C2a_3 = 3C2A_3) Secondary.siRNA.Screen.Data\(siRNA_Target <- recode(Secondary.siRNA.Screen.Data\)siRNA_Target, PIK3C2a_4 = 3C2A_4) cell.numbers.secondary.screen <- count(Secondary.siRNA.Screen.Data, siRNA_Target) count(Secondary.siRNA.Screen.Data, siRNA_Target) r mean.cell.numbers.secondary.screen <- summarise(cell.numbers.secondary.screen, avg = mean(n), SD = sd(n), sum(n)) summarise(cell.numbers.secondary.screen, avg = mean(n), SD = sd(n), sum(n))
rr svglite(paste(panel_output_path, Secondary_siRNA_Screen_boxplots_filename, .svg, sep = \), width=7, height=5) ggplot(Secondary.siRNA.Screen.Data, aes(x=siRNA_Target, y=Ratio_A.CMAC.to.T.CMAC_mean_intensity.z.score, color = factor(siRNA_Treatment_Type))) + geom_boxplot(notch = TRUE, outlier.colour = ) + theme(axis.text.x = element_text(angle=90)) + ylim(-3, 6) + scale_x_discrete(limits=c(,,_3Ctrl_Pool, -TARGETplus Non-targeting Control, Non-targeting Control, , , 5, _1, _2, _3, _4, 4K2A_1, 4K2A_2, 4K2A_3, 4K2A_4, 4KA_1, 4KA_2, 4KA_3, 4KA_4, 5K1B_1, 5K1B_2, 5K1B_3, 5K1B_4, 5K1C_1, 5K1C_2, 5K1C_3, 5K1C_4, _1, _2, _3, _4, 3C2A_1, 3C2A_2, 3C2A_3, 3C2A_4)) dev.off()
null device
1
rr ggplot(Secondary.siRNA.Screen.Data, aes(x=siRNA_Target, y=Ratio_A.CMAC.to.T.CMAC_mean_intensity.z.score, color = factor(siRNA_Treatment_Type))) + geom_boxplot(notch = TRUE, outlier.colour = ) + theme(axis.text.x = element_text(angle=90)) + ylim(-3, 6) + scale_x_discrete(limits=c(,,_3Ctrl_Pool, -TARGETplus Non-targeting Control, Non-targeting Control, , , 5, _1, _2, _3, _4, 4K2A_1, 4K2A_2, 4K2A_3, 4K2A_4, 4KA_1, 4KA_2, 4KA_3, 4KA_4, 5K1B_1, 5K1B_2, 5K1B_3, 5K1B_4, 5K1C_1, 5K1C_2, 5K1C_3, 5K1C_4, _1, _2, _3, _4, 3C2A_1, 3C2A_2, 3C2A_3, 3C2A_4))
rr siRNA.conditions.for.U.test <- subset(Secondary.siRNA.Screen.Data, siRNA_Target != & siRNA_Target != & siRNA_Target != _3Ctrl_Pool & siRNA_Target != Non-targeting Control) p_values <- c() for (i in 1:length(unique(siRNA.conditions.for.U.test\(siRNA_Target))){ for (j in 1:length(unique(siRNA.conditions.for.U.test\)siRNA_Target))){ siRNA_Target_1 = subset(siRNA.conditions.for.U.test, siRNA.conditions.for.U.test\(siRNA_Target == unique(siRNA.conditions.for.U.test\)siRNA_Target)[i]) siRNA_Target_1\(siRNA_Target_Number = as.factor(i) siRNA_Target_2 = subset(siRNA.conditions.for.U.test, siRNA.conditions.for.U.test\)siRNA_Target == unique(siRNA.conditions.for.U.test\(siRNA_Target)[j]) siRNA_Target_2\)siRNA_Target_Number = as.factor(j+100) siRNA_Target_Pair = rbind(siRNA_Target_1, siRNA_Target_2) siRNA_Target_Pair_ttest = wilcox.test(Ratio_A.CMAC.to.T.CMAC_mean_intensity.z.score ~ siRNA_Target_Number, data = siRNA_Target_Pair)\(p.value p_values = rbind(p_values, siRNA_Target_Pair_ttest) } } p_values <- as.vector(p_values) p_value_matrix <- matrix(p_values,nrow = length(unique(siRNA.conditions.for.U.test\)siRNA_Target)),ncol = length(unique(siRNA.conditions.for.U.test\(siRNA_Target))) columns <- as.vector(unique(siRNA.conditions.for.U.test\)siRNA_Target)) rows <- as.vector(unique(siRNA.conditions.for.U.test\(siRNA_Target)) p_value_dataframe <- as.data.frame(p_value_matrix, row.names = rows) names(p_value_dataframe) <- columns # print(p_value_dataframe) # Introduce Holm (also called Holm-Bonferroni) p-value correction (used because it gives same stringency against false positives (type 1 errors) with lower probability of false negatives (type 2 errors))) corrected_p_values <- p.adjust(p_values, method = \holm\, n = length(p_values)) corrected_p_value_matrix9 <- matrix(corrected_p_values,nrow = length(unique(siRNA.conditions.for.U.test\)siRNA_Target)),ncol = length(unique(siRNA.conditions.for.U.test\(siRNA_Target))) columns <- as.vector(unique(siRNA.conditions.for.U.test\)siRNA_Target)) rows <- as.vector(unique(siRNA.conditions.for.U.test$siRNA_Target)) corrected_p_value_dataframe <- as.data.frame(corrected_p_value_matrix9, row.names = rows) names(corrected_p_value_dataframe) <- columns print(corrected_p_value_dataframe)
rr heatmap.2(log(corrected_p_value_matrix9), dendrogram = ‘none’, Rowv = FALSE, Colv = FALSE, labRow = rows, labCol = columns, srtCol = 45, cexRow = 1, cexCol = 1, margins = c(7, 9), symm = TRUE, revC = FALSE, breaks = c(log(1e-300), log(1e-50), log(1e-10), log(1e-6), log(0.001), log(0.01), log(0.05), log(1)), rowsep = 1:nrow(corrected_p_value_matrix9), colsep = 1:ncol(corrected_p_value_matrix9), sepcolor = ‘darkgrey’, sepwidth = c(0.02,0.02), trace = ‘none’, col = c(, , , , , , ), denscol = NULL, keysize = 1.5, key.title = NA)